Local Moment Formation and Kondo Effect in Defective Graphene 
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We study the local moment formation and the Kondo effect at single-atom vacancies in Graphene. 
We develop a model accounting for the vacancy reconstruction as well as non-planarity effects in- 
duced by strain and/or temperature. Thus, we find that the dangling a orbital localized at the 
vacancy is allowed to strongly hybridize with the 7r-band since the scattering with the vacancy turns 
the hybridization into singular function of the energy (~ [|e| In 2 e/D]' 1 , D ~ the bandwidth). This 
leads to several new types of impurity phases, which control the magnitude of the vacancy magnetic 
moment and the possibility of Kondo effect depending on the strength of the local Coulomb inter- 
actions, the Hund's rule coupling, the doping level, and the degree of particle-symmetry breaking. 



PACS numbers: 

Experimental evidence of magnetic local moments 
(MLM) in irradiated graphene at very low temperatures 
has been recently reported in [TJ [5] . In another experi- 
ment using irradiated samples [3J, the Kondo effect with 
surprisingly high Kondo temperatures (Tk ~ 10-100 K, 
obtained by fitting resistivity data to conventional spin- 
1/2 Kondo laws [4]) was observed. A 'priori, the latter 
observation appears to be at odds with the observations 
of [TJ |5] ■ For instance, at the temperatures of the Manch- 
ester group experiment (T = 2K) [T], the measured spin- 




FIG. 1: Two level Anderson model describing a reconstructed 
single-atom vacancy in rippled graphene; d is an orbital that is 
mainly localized at a dangling a orbital of the apical atom; ttq 
is the vacancy-induced localized state. Rippling implies that 
the d-orbital can hybridize with the n conduction band via 
sp 2 -sp 3 hybridization, V. The hybridization, A(e), is energy 
(e) dependent and enhanced by the strong scattering potential 
at the vacancy. U a7V and Jh describe the Coulomb interaction 
between the two localized levels. See Eq. p). 



1 /2 MLM of the vacancies created by irradiation should 
be quenched by the Kondo effect observed in [3J. Fur- 
thermore, the authors of [3J reported a Tk with a rather 
weak dependence on gate voltage. However, theories of 
the Kondo effect in graphene [5H5] predict a critical cou- 
pling below which MLM are stable at the Dirac point 
(DP) whilst Tk remains small around the DP [318]. 

The puzzle described above calls for re-visiting the 
magnetic properties of isolated single atom vacancies in 
graphene. Here we show that earlier theories [6-9 J must 
be modified to account for the strong scattering poten- 
tial of the vacancy. Indeed, the later is responsible for 
a dramatic change of the energy dependence of the hy- 
bridization with the 7r band of a level localized at a recon- 
structed vacancy [T2"HT4"] . Thus, the energy e dependence 
of the hybridization is changed to ~ [|e| ln 2 (e/D)] _1 near 
the DP (cf. Fig. [T]), which is different from the pseu- 
dogap behavior ~ |e| r (r > 0) considered earlier |6r[9]. 
Therefore, it becomes possible for vacancy MLM to be 
quenched at T K ~ 1 - 100 K (cf. Fig|J. 

The possibility of defect-induced magnetism in 
graphene has been widely discussed in the literature (see 
e.g. [lTl - [T3l [15], and references therein), mainly using 
various types of static mean-field theories. However, a 
recent study [15] using dynamical mean-field theory con- 
cluded that vacancies should allow for the formation of 
MLM that are ferromagnetically coupled to the conduc- 
tion band of graphene. Within Hartree-Fock, this hap- 
pens because of the strong local correlations [16] and the 
large enhancement of the local density of states in the 
neighborhood of the vacancy [T7]- However, the anal- 
ysis of Refs. [TTJ [TS] focuses only on the contribution 
to the vacancy magnetic moment from 7r-band electrons. 
Ab initio calculations using the generalized gradient ap- 
proximation (GGA) to density functional theory (DFT) 
[T2l IT5] , can account for the contribution of the dangling 
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a orbitals and are also able to describe the possible re- 
constructions of the vacancy [ID] . However, such studies 
predict different values for the vacancy MLM [L2l - fl4] . 



Indeed, provided the dangling a orbitals are not passi- 
vated, the vacancy reconstructs [TQl [12j H3J [18] following 
a Jahn- Teller distortion. Thus, a carbon pentagon with 
a strong bond resulting from two of the dangling a or- 
bitals (cf. Fig. [I]) is formed. The remaining dangling a 
orbital is mainly localized at the apical atom (opposite 
to the reconstructed bond) , and appears as a level at ~ 1 
eV [H [13] below the DP (cf. Fig.Q. Double occupation 
of this orbital is strongly suppressed by a rather strong 
Coulomb repulsion U ~ 10 eV [TBJ, which leads to the 
formation of a MLM. However, in flat graphene, symme- 
try forbids the hybridization of a and 7r orbitals, mean- 
ing that the electron in the dangling a orbital cannot 
hybridize with 7r band and thus exchange with conduc- 
tion electrons. In absence of other contributions to the 
MLM [12 , this fact is hard to reconcile with the observed 
vacancy effects on charge [3J and spin transport 0. 

The above picture changes if we recall that graphene is 
a membrane that becomes easily rippled under the effect 
of strain and temperature [19]. Indeed, DFT-GGA cal- 
culations [18] indicate that, under rather small (ps 1%) 
isotropic compression, regions containing reconstructed 
vacancies do ripple [20 . Rippling is achieved by means of 
sp 2 -sp 3 hybridization of the carbon bonds [3TJ [55] , which 
admixes the a with p 2 orbitals at a given atom, resulting 
in new a orbitals pointing along directions that are no 
longer perpendicular to the p z orbital. Thus, in a typi- 
cal sp 2 -sp 3 hybridized configuration, a carbon atom can 
stick out of the graphene plane. If the atom contains a 
dangling a orbital, as in the case of the apical atom of 
a reconstructed vacancy, the a orbital will be able hy- 
bridize with the 7r band. In fact, the hybridization en- 
ergy, V a7T , can be substantial already for small deviations 

from planarity [HI HI]: V av = A^J±=^-(e s - e p ), where 

sin# = A/\f A 2 + 2, 9 being the angle subtended by the 
dangling a orbital atom and the graphene plane, and 
e s — e p = — 8.31eV is the atomic s-p splitting. Thus, in 
rippled graphene tt electrons are allowed to hop on and 
off the dangling a orbital, thus leading to an effective 
(anti-ferromagnetic, AF) exchange with the conduction 
band. However, this also means that the MLM of the 
dangling er orbital may be quenched by the Kondo effect. 
Indeed, this possibility is not negligible, since the strong 
scattering potential of the nearby vacancy substantially 
modifies the local density of states [17] . 

Following the previous discussion, we introduce a 
model that contains all the ingredients described 
above 1231 and whose Hamiltonian can be written as 



H = H d + H hyb + H V+7T , where [3D] 

H d = ( d n d + Un d ^-n di + H dn , 
Fhyb = V [d a b a + bU a ] 

H v+7! = e Q al a% - t ^ a\ a b^ + 



(1) 
(2) 
(3) 



The term H d above describes the dangling a orbital 
with e d ~ —0.75 eV measured from the DP [ID] , and 
U w 10 eV P3J; n dt(l) = dL^tfM ; H dlT (see below) 
describes the Coulomb interaction between the localized 
dangling a orbital (d-orbital, in what follows) and the 7r 
band electrons below. Hh y b is the hybridization with the 
7r-band, where V = M°W^ (M0)| 2 =s 0.7 [13 [23]). 
In the last term , the limit £o —¥ +oo must be taken in 
order to describe a single-atom vacancy on the A sublat- 
tice at R = (0, 0). The apical atom corresponds to the B 
sublattice atom within the same lattice. The fermion op- 
erators af, bf (a ai) b^A destroy (create) 7r-band electrons 
of spin a on sites belonging to A and B sublattices. 

To begin with, we retain only the nearest neighbor 
hopping (i ~ 2.8 eV [5j) and consider the chemical 
potential to be at the DP (i.e. fi = 0). The ef- 
fect of a nonzero /z and the next-nearest neighbor hop- 
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af + b 



t b a 

la 3 



will be discussed at the 



ping -t' J2 

end. Thus, the spectrum of ff^+v contains a localized 
level (a zero-mode, ZM) pinned at the DP, which is or- 
thogonal to a continuum of waves scattering off the va- 
cancy |17j . The ZM wavefunction is not square normal- 
izable since |y>o(R) 1 2 ~ ]Kp' ^" u ^ ^ s over l a P with the 
d-orbital is the largest, which means that the dominant 
contribution to H dlT is from the Coulomb interaction be- 
tween the electrons in the ZM and the <i-orbital [23], i.e. 
H dlT = U an n d n Q - J H Sd ■ So, where n d = + n di , 



no 



3 7Tq , where 7Tq (Tr a0 ) destroys (creates) electrons 



P 



in the ZM, and S d = 4 (f ) J dP, and S = tt 0q (f ) fi „ , 
where cr = {a x 1 a v , a z ) are the Pauli matrices. Jh > is 
the Hund's coupling and U aT > is the Coulomb repul- 
sion between electrons in the ZM and the d orbitals. The 
relative value of these parameters determines whether the 
vacancy spin is a doublet or a triplet. However, the ex- 
tended nature of the ZM makes it difficult to obtain accu- 
rate estimates for Jjj and U alT (some estimates, obtained 
on dimensional grounds, arc provided in [23 ). Thus, in 
absence of accurate input, below we explore the possible 
magnetic phases for both the doublet and the triplet. 

Treating electrons in scattering states as non- 
interacting, allows us to integrate them out exactly. 
Their influence on the <i-level can be encapsulated in 
a self-energy function [I4|, whose imaginary part is 
7v\V\ 2 p B (e), where p B (e) = p B (e ^ 0), being p B (e) the 
local density of states at the apical atom [53]. The latter 
exhibits a singular behavior p B (e) ~ m ln 2 1 | e / £) | (D ~ t) 
around the DP [17] . We can mimic the effect of the self- 
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energy by representing the it band as a electron reser- 
voir with energy-dependent tunneling into the d-level [4] , 
leading to 



r+D 

H = H v + H hYh + det^ a (e) </> Q (e) , 

J-D 

H v = edUd + Und^ridi + e D n + U,y^n d n - J H Sd ■ So, 
H hyb = V j +D de [t(e)dlr (e) + t* (e)V>£ (e)d a ] , (4) 

J-D 

where e D = 0, A(e) = \t{e)\ 2 = ps(e), and D = D - i is 
a high-energy cut-off of order of the width of the ir band. 
Since (for small angle 0) /S.{e d )/\£d\ <C 1, charge fluctu- 
ations in the dangling er orbital are small at low tem- 
peratures. However, the charge state of the ZM is deter- 
mined by the magnitude of the parameter G = Jh/^Uo-^: 
(rid) ~ (doublet) for G -C 1 , and (rid) ~ 1 (triplet) for 
G ^> 1. In both cases, the charge fluctuations can be inte- 
grated out using a Schrieffer- Wolff transformation, which 
leads to a Kondo model: 



H K = H C + J K s ■ S v + Voflofg, 



(5) 



where s = 4„ (f /o Q = I det(e) ^(e) and S v 

is the total spin operator of the vacancy: Sy = S^, for 
the doublet (G <C 1) and Sy = S d + So, for the triplet 
(G» 1). The expressions for Jk and the potential Vq 
depend on the total spin of the vacancy [23] . 

However, for G sa 1, the charge fluctuations in the 
ZM level are not negligible and the description in terms 
of Eq. ^ breaks down. In such a case, we must deal 
with the full Anderson model, Eq. Q. Having G « 1 
requires fine tuning as it corresponds to a quantum phase 
transition between the doublet and the triplet, which will 



be analyzed elsewhere 



Here, we assume that U„~ 1 1 - 



G\ is large enough so that using Eq. ^ with a cut-off 
Do < minder |, c/ CT7r |l — G\} is sufficient. 

Next, we carry out a scaling analysis of ^ [4l IS] 125], 
The perturbative scaling equations for jk = Jk /Dq and 
v = ^o/A) read: 



dj 



K 



dhiD 



= r(D)j K -f K + 0(f K ), 



dv 



d\nD 



r(D)v, (6) 



where r(D) = (din A(e)/dlne) £= £> is assumed to be a 
slowly varying function of D, which is the case since near 
the DP A(e) - [|e| In^/D)}' 1 and r(D) ~ -1 with loga- 
rithmic accuracy. Thus, as D is reduced, jx(D) and v (D) 
grow and the vacancy enters the strong coupling regime. 
For jx{D), a similar conclusion was reached in [25] based 
on an analysis of the orthogonality catastrophe in defec- 
tive graphene. In the strong coupling regimem where 
Jk{D*) ~ 1 (or v(D*) ~ 1), ^ break down. However, 
the Kondo temperature can be related to the crossover 
scale D* where jk{D*) ~ 1 (provided jk ~ 1 before 
v). The solution of @ yields T K \n(T K /D ) ~ J K 



Vacancy Spin 


jk ~ 1 


V ~ 1 


Doublet (Sv = |) 


ssc 


ALM 


Triplet (Sv = 1) 


ussc 


ALM 



TABLE I: Stable phases for the Kondo model of Eq. {5|. SSC 
stands for singlet strong coupling , ALM for asymmetric local 
moment, and USSC for underscreened strong coupling; jk ~ 
1 and v ~ 1 indicate which coupling (either jk oc Jk or 
becomes relevant first. 



for Jx <C -Do -•'> • Note that, this means a dramatic 
enhancement of Tk compared to the vanishing (expo- 
nentially small) Tjc of the pseudo-gap [6] [7j [9] [25] (flat- 
band [4]) Kondo models. 

In the strong coupling regime, we have to rely on other 
methods different from Eqs. In this regard, there are 
remarkablesimilarities of the present Kondo model with 
one studied by Bulla and Vojta [37] using the numeri- 
cal renormalization group (NRG). These authors consid- 
ered a spin-I/2 Kondo model where A(e) ~ |e| r , with 
— 1 < r < 0. With logarithmic accuracy, we can bor- 
row their results for r \ — 1 and Jk > 0. Thus, two 
stable strong coupling fixed-points (FPs) are known to 
exist: a particle-hole asymmetric local moment (ALM) 
FP, at which the vacancy magnetic moment is decou- 
ples from the band, and a singlet strong-coupling fixed 
point (SSC) FP, at which the single electron in dangling 
a orbital forms a singlet with the conduction band elec- 
trons, and vacancy magnetic moment is thus quenched. 
The flow between ALM and SSC FPs is controlled by a 
critical point |27j . 

For the triplet case (G> I), to the best of our knowl- 
edge, there are no NRG results available. However, in 
this case we can tentatively assume that, as jx flows 
to strong coupling, half of the vacancy spin- 1 will be 
quenched. Applying an argument due to Nozieres and 
Blandin [55], the residual spin- 1/2 couples ferromagneti- 
cally (FM) to the conduction band, that is, for T <C Tk , 
it is described by Hsc = H c + J'si • S' v where S' v is the 
residual vacancy spin, and J 1 < 0. In this strong-coupling 
coupling picture, one electron is captured at the first site 
(site in our notation above) of the Wilson chain [J] rep- 
resentation of H c (cf. Eq. thus projecting out this 
site. The residual spin- 1/2 couples to the next site of 
the Wilson chain, whose spin is described by the oper- 
ator sj.. The hybridization function for this site can be 
obtained by tridiagonalization, which allows to relate the 
local Green's functions at sites and 1 with site pro- 
jected out (see g|, chap. 4): gi (e) = - |y|2 ^ (e) , 

with V cx V. Hence A x (e) cx -Im gi (e) - \e\/ In 2 \e/D\. 
Hence, with logarithmic accuracy, r(D) ~ 1 in Eq. |6j), 
which, together with j_r-(0) ~ J' < 0, implies that J' 
renormalizes to zero. Thus, Hsc describes a stable un- 
derscreened Kondo FP. Moreover, as in the double case, 
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FIG. 2: For a vacancy with S v = §, large N mean-field 
Kondo temperature Tk vs. Jx/Fd (Jk = 2jCj 2 /|e<j and 
r d = tt| V| 2 A(e d ),(ed = -0.75 eV, t = 2.97 eV and (7 -> +oo). 
Panel (a): /i = and t' — (continuous), /i = 7.4 meV and 
£' = (dashed), and fj, = and t' = 30 meV. Panel (b): 
i' = 59 meV and \i = —30 meV (continuous), fj, = —22 meV 
(dotted), fi = —37 meV (dot-dashed). The fi, t' = result is 
also plotted for reference (dashed). Note that for Jk/Tci ~ 1 
the vacancy enters the mixed-valence regime. 



for large particle- hole asymmetry, Vq, it is also reasonable 
to expect the existence of a strong-coupling fixed similar 
to the ALM FP, and a critical point in the V-Jk plane 
that controls the flow between the two FPs. The possible 
phases of ([5| are summarized in Table [IJ 

Finally, let us consider the effect and the next neigh- 
bor hopping, t' , and a small doping /i / 0. At first 
sight, both perturbations introduce particle-hole asym- 
metry |17j . which, for fi,t' <C t, can be regarded as ad- 
ditional contributions to the bare value of Vo in Eq. ( pjj) . 
Thus, the initial conditions for the scaling equations Ml) 
are modified from their Schricffer- Wolff values , which 
favors the flow towards the ALM FP. For /i 7^ 0, the flow 
towards the ALM FP can be seen as the first stage of 
a two-stage flow [5], whose final stage occurs at much 
lower temperatures and may be a crossover from the 
ALM to a more conventional Kondo-singlet FP (see be- 
low). The RG flow for the triplet case will be analyzed 
elsewhere [24] . 

To test the above ideas and obtain some quantitative 
estimates for Tjc , we have used a slave-boson mean-field 
approach to Anderson model in Eq. Q in the U — > +00 
limit HI2S]. In the mean-field theory [41155]. Tk is de- 
termined from the equation: 



T, 



K 



g(Wn) e to n 0+ 



1 

Jk'' 



(7) 



where w„ = 2n(n + \)T K , J K = 2\V\ 2 /\e d \ is the 
Kondo coupling for U — > +00 [23], and g(iio n ) — 
f deA(e)/(iuj n — e + fx) is the local Green's function at 
the apical atom. The results of a numerical solution of 
this equation are displayed in Fig. [2] From Fig. [2ja), it 



can be seen that turning on a finite \i or t' strongly sup- 
presses Tk- For finite t', Tk appears to vanish at small 
Jk <C Tk{Jk,H = 0, t' 7^ 0), which we interpret as the 
ALM FP being favored over the SSC FP. However, the 
effect of /j, 7^ seems to be slightly different because the 
curve for Tk changes its concavity, suggesting a finite but 
small Tk at small Jk (finite-size effects prevent us from 
obtaining accurate results for Tk below ~ 1 K). Fig[2jb) 
displays the values of T K for t'/t = 0.02 (i.e. t' = 59 
meV) and several values of /i. We find that when t has 
opposite sign to fi, their particle-hole symmetry-breaking 
contributions appear to partially cancel each other, lead- 
ing to a less pronounced suppression of Tk (relative to 
the results in Fig.[2]ja)). We believe this is because a re- 
duction of the bare value of Vo in Eq. ^ favors the SSC 
over the ALM FP for large enough Jk- 

In conclusion, we have studied a model for a recon- 
structed single-atom vacancy in graphene, finding that 
the strong vacancy scattering potential can dramatically 
affect its magnetic properties. We find that the, if 
particle-hole symmetry breaking is not too strong, Kondo 
effect can occur for small doping levels provided the dan- 
gling a orbital is not passivated. However, for sufficiently 
large particle-hole symmetry breaking, a local moment 
that becomes increasingly decoupled from the n band as 
temperatures decreases appears. Further consequences 
of our theory for existing and future experiments will be 
presented elsewhere |24) . 

We acknowledge illuminating discussions with V. 
Pereira, I. Grigorieva, D. Sanchez-Portal, and E. J. G. 
Santos. MAC and FG acknowlege the hospitality of 
KITP, Santa Barbara. This research has been supported 
by ERC, grant 290846 and MICINN, Spain, (FIS2008- 
00124, FIS2011-23713, CONSOLIDER CSD2007- 00010) 
(FG) and FIS2010-19609-C02-02 (MAC), and partially 
supported by the National Science Foundation under 
Grant No. NSF PHY11-25915. AI acknowledges fi- 
nancial support from CONICET (PIP 0662), ANPCyT 
(PICT 2010- 1907) and UNLP (PID X497), Argentina. 



[1] R. R. Nair, M. Sepioni, I.-L. Tsai, O. Lehtinen, J. 

Keinonen, A. V. Krasheninnikov, T. Thomson, A. K. 

Geim, I. V. Grigorieva, Nat. Phys. 8, 199 (2012). 
[2] K. M. McCreary, A. G. Swartz, W. Han, J. Fabian, 

R. K. Kawakami arxiv: 12062628 (2012). 
[3] J.-H. Chen, L. Li, W. G. Cullen, E. D. Williams, and 

Michael S. Fuhrer, Nat. Phys. 7, 535 (2010). 
[4] A. C. Hewson, The Kondo Effect to Heavy Fermions 

(Cambridge University Press, Cambridge, 1993). 
[5] A. H. Castro-Neto, F. Guinea, N. M. R. Peres, K. S. 

Novoselov, and A. Geim, Rev. Mod. Phys. 81, 109 (2009). 
[6] F. Withoff and E. Fradkin, Phys. Rev. Lett. 64 1863 

(1990). 

[7] K. Sengupta and G. Baskaran, Phys. Rev. B 77, 045417 
(2007). 



5 



[8] M. Vojta, L. Fritz, R. Bulla, Europhys. Lett. 90, 27006 
(2010)'. 

[9] B. Uchoa, T. G. Rappoport, and A. H. Castro Neto, 
Phys. Rev. Lett. 106, 016801 (2011). 
[10] A. Zunger and R. Englman, Phys. Rev. B 17 (1978); H. 

Amara et al, Phys. Rev. B 76, 115432 (2007). 
[11] M. P. Lopez-Sancho, F. de Juan, M. A. H. Vozmediano, 
Phys. Rev. B 79, 075413 (2009) 



[12] J. J. Palacios and F. Yndurain, arXiv:1203.6485 (2012). 
[13] N. Nanda et al, arxiv: |arXiv:1105.1129| (2011), to appear 

in New J. of Phys. 
[14] O. V. Yazyev and Lothar Helm, Phys. Rev. B 75 125408 
(2007). 

P. Haase, S. Fuchs, T. Pruschke, H. Ochoa, and F. 
Guinea, Physical Review B 83, 241408 (2011). 
J. A. Verges, E. SanFabian, G. Chiappe, G. Chiappe, 
and E. Louis, Phys. Rev. B 81, 085120 (2010); T. O. 
Wehling, E. Sasioglu, C. Friedrich, A. I. Lichtenstein, M. 
I. Katsnelson, and S. Bhigel, Phys. Rev. Lett. 106, 236805 
(2011); 

[17] V. M. Pereira et al, Phys. Rev. Lett. 96, 036801 (2006); 

N. M. R. Peres, S. W. Tsai, J. E. Santos, and R. M. 

Ribeiro, Phys. Rev. B 79, 155442 (2009). 
[18] E. J. G. Santos, S. Riikonen, D. Sanchez-Portal, A. 



[15] 
[16] 



Ayuela, Journal of Physical Chemistry 116, 7602 (2011) 
E. J. G. dos Santos, PhD Thesis, The University of the 
Basque Country, Spain (2011). 

[19] J. C. Meyer et al, Nature 446, 60 (2007); A. Fasolino, 
J. H. Los, and M. I. Katsnelson, 6, 858 (2007); W. Bao 
et al, Nature Nanotech. 4 562 (2009). 

[20] Indeed, curbed carbon allotropes such like Ceo also con- 
tain pentagons. 

[21] L. Pauling, The Nature of the Chemical Bond (Cornell 
University Press, Ithaca, NY, 1972). 

[22] A. H. Castro-Neto and F. Guinea, Phys. Rev. Lett. 103, 
026804 (2009). 

[23] See supplementary material, available from this URL. 
[24] M. A. Cazalilla, A. Iucci, A. H. Castro Neto, and F. 

Guinea , in preparation. 
[25] L. Dell'Anna, J. of Stat. Mech, P01007 (2010). 
[26] M. Hentschel and F. Guinea, Phys. Rev. B 76, 115407 

(2007). 

[27] M. Vojta and R. Bulla, Eur. Phys. J. B 28, 283 (2002). 
[28] P. Nozieres and A. Blandin, J. Physique 41, 193 (1980). 
[29] D. News and N. Read, Adv. Phys. 36, 799 (1987). 
[30] The convention of summing over repeated Greek indices 
(a =f, i) will be used throughout. 



large Supplementary Material 



DETAILS OF THE MODEL OF THE VACANCY 



We next describe the derivation of the model that is analyzed in the main text. In establishing a minimal model for 
a reconstructed single-atom vacancy in graphene, we shall rely on existing first-principle calculations (e.g. Refs.[THI|- 
The results of those calculations can be summarized as follows: 

1. Removing one carbon atom (i.e. creating a vacancy) leaves three dangling er-bonds behind. The latter manifest 
themselves as three impurity levels splitting off the sigma band. 

2. From a tight-binding perspective, the vacancy cr-levels can be regarded as essentially the dangling a bonds of 
the atoms around the vacancy plus some admixture of a states away from the vacancy [2] . 

3. In an unrelaxed lattice, where all carbon atoms around the vacancy site remain at the same distance as the 
perfect graphene lattice, the crystal field effects, which correspond to the nearest neighbor hopping between 
cr-orbitals in the original lattice) split the tree dangling-bond levels into one lower energy level (the symmetric 
combination) and two degenerate excited levels (cf . Fig. [1| . 

4. The situation described in the previous point is energetically unstable against a Jahn- Teller distortion, which 
breaks the D$h symmetry of the unrelaxed lattice and draws two of the three carbon atoms surrounding the 
vacancy closer to each other (cf . Fig. [I] in the main text and and Fig. |3| . 

5. Further ah initio calculations [5] suggest that, in the presence of isotropic compression, the three atoms around 
the impurity are no longer coplanar. This is possible if e.g. the apical atom (which contains a dangling a- 
bond [2 ) is lifted slightly above (or below) the plane of the two atoms opposite to it (i.e. those that form a 
reconstructed a bond as a result of the Jahn- Teller distortion), see Fig. [T] in the main text. This non-coplanarity 
is made it possible by sp 2 -sp 3 hybridization of the atoms around the vacancy. In the model below, we assume 
that such hybridization is particularly important in the case of the apical atom (cf. Fig. [l]in the main text). 
This allows the electrons to hop in and out of the dangling a orbital at the apical site. This effect can be also 
achieved by rippling that occurs at any finite temperature. 
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Given the above considerations, we propose the following model of the vacancy: [TU] 



H 



7T+V 



E Htf 



e E a aO a aO, 



(id), 
3 



ffhyb = E [ 6 lorf? + &L>#] + V m E 6 L> [&L>«? + fo lo4 + 4l&0 + 4&2 



+ #V-int) 



(1) 

(2) 

(3) 
(4) 



where T > T' is the hopping between the vacancy tr-orbitals and V^o- describes the sp 2 -sp 3 hybridization [5] of the 
apical atom (see Fig. [3] for details about the labeling of the different orbitals): 



V m = A 



I- A 2 



(5) 



where sin# = Aj\J A 2 + 2, 9 being the angle between the a orbitals at the apical atom site (cf. Fig. |3j) and the 
Graphene plane, and e s — e p — — 8.31eV is the atomic s-p splitting. The above model neglects the local modifications 
of the hopping parameters around the vacancy site, both resulting from the Jahn- Teller distortion and from the sp 2 - 
sp 3 hybridization caused by an isotropic compression. In Eq. Q Hd n describes the Coulomb interaction between the 
a and the 7r orbitals, which is particularly strong for the 7r and a orbitals localized around the vacancy (see Sec. |). 

Furthermore, In Eq. ([I]), the term Hy-int describes the interactions between the electrons localized at the vacancy-a 
orbitals, which we will specify below. In the above equation, the limit e ~ ^ +°° is implicitly assumed in order to 
remove the site at R = (0, 0) of the A sublattice. 

As discussed in point [4] above, two of the dangling bonds bind when the lattice undergoes a Jahn- Teller distortion. 
If we diagonalize the subspace spanned by the vacancy a orbitals by means of the transformation: 




(u n ,u 12 ,u 13 ) 
, with r 



l + %/l+8r 2 
1r 



( U\\ Ui2 Uiz 

U 21 U 22 U 2 3 

\ U 3 x U 32 "33 

,1)1), U 2 = (U21,U 2 2,U 23 ) = J\f 2 { 




(6) 



l-Vl+8r 2 
2r 



1, 1), u 3 = («31,«32,«33) 



where U! 

77g(0, 1, —1), with r = T'/T and Mi the normalization constants. The energies of the levels resulting from the previous 
diagonalization are given by first term in 



Hy — E £i ^L^i* + Hv-int> 



(7) 



a, i—1 



< B 

s 2 # Si 

A° R= (0,0) ' 

• v 

d 3 d 2 



FIG. 3: Labeling of the relevant orbitals around the vacancy. We make the approximation that the vacancy <r-orbitals correspond 
essentially to the dangling bonds around the vacancy site (empty circle). See Eq. Q for details. 
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where ei = e CT + §T (\/l + 8r 2 - l), e 2 = e CT - §T (\/l + 8r 2 + l), e 3 = e CT + T (According to the fitting to DFT- 
LDA calculations carried out by Nanda et al. in Ref. [2j T — 1.85 eV and T' = 1.30 eV). Thus, the level number 
3 corresponds to an excited state, which lies above the Dirac point, [2] and the level 2 corresponds to a deep level 
describing the reconstructed bond at the basis of the pentagon opposite to the apical atom (cf. Fig. [I] in the main 
text and[3|. 

The remaining electron from the vacancy a orbitals occupies the orbital with energy ei, and in flat graphene it 
cannot not hybridize with the 7r-band although it has a very large amplitude at the apical atom (see expression for 
Ui above). However, in rippled graphene (as in the presence of isotropic compression) there is hybridization, through 
the term -f/hyb m Eq. Q. In the following, we shall investigate whether such hybridization could lead to a Kondo 
effect. Thus, we shall focus on the orbital 1 and, to lighten the notation, henceforth we shall drop the orbital index 
(the other two orbitals, 2 and 3, will be treated as inert). Furthermore, we shall neglect the last term in Eq. Q, that 
is, the sp 2 -sp 3 hybridization with the s\ and s 2 er-orbitals of the apical atom (cf. Fig. p3J). Such terms are expected to 
lead to a renormalization of the hybridization of the vacancy level 1 with the 7r band. These considerations lead to 
the following Anderson model: 

a a 

where H^ +v is given by Eq. Q with e ->■ +oo; e d = e x = e CT + \T (a/I + 8r 2 - l) ~ -0.75 eV, [2j |4| and V — 
\ipd(0)\ 2 V a7: (6), where |(^d(0)| 2 w 0.7, according to the above tight-binding model, and U « 10 eV [7]. 

The novelty of ([8| lies in the fact that it is a single-impurity Anderson model in the presence of a strong scattering 
potential that induces a resonance at the Dirac point. Note that in the above model we have neglected any interactions 
between the electrons in the 7r-band. 



e F 




Vacancy Crystal Jahn-Teller 
a levels field distortion 

FIG. 4: Vacancy a-levels and evolution when different perturbations are taken into account. For left to right: i) Crystal field 
splitting (i.e. nearest neighbor hopping T — T' and ii) Jahn-Teller distortion with T > T' . The lowest level, which has equal 
weight in the three atoms surrounding the vacancy has lowest energy and it is assumed to be doubly occupied, thus forming a 
(stretched) a bound at the basis of a pentagon (cf. figure [l] in the main text). The intermediate level is below the Fermi energy 
€f (assumed to be at the Dirac point) and we shall assume it to be singly occupied as U ~ 10 eV [7] 



TWO LEVEL ANDERSON MODEL 

Let us return to the above Hamiltonian and consider first the scattering problem of the tt band electrons with the 
vacancy. Let us assume that we have diagonalized the Hamiltonian ^ and therefore we can write: 

H-n+V = £07ri 71"0 + ^2 e pn ^pan^pn^ ( 9 ) 
a,p,n 

where the index p = c,v describes the valence and conduction band scattering channels and a =t, I- The other 
quantum numbers have been collectively denoted by n. We have explicitly separated the zero mode, whose energy 



(8) 
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coincides with the Dirac point energy (eq = for a particle-hole symmetric case where we neglect next nearest neighbor 
hopping) , and it described by TV a , Tr a . The original lattice operators can be expressed in terms of the eigenmodes of 
H^+v as follows: 

a? = ]>> Apn (R 4 )^ n , (10) 

p . n 

bf = po(Ri)7tf + J2 fBpni^) (11) 
p,n 

In the last expression above, we have used the fact that the zero-mode wavefunction is (within the tight-biding 
nearest-neighbor hopping approximation) is fully contained within the majority B sublattice (see e.g. Ref. EH & n d 
references therein). 

In this eigenmode basis, the minimal Anderson Hamiltonian, Eq. ([8]), introduced in previous section reads: 
H AM = H« +v + ed^dld" + Ud\d\d^ + V £ [<A)(0)4< + ^(0)4o^ Q 

a a 
+ V [fB P n(0) 4<n + V*B P n(0X P nd a ] + H d7r (12) 

(13) 

Note that the zero-mode wavefunction is not normalizable since (neglecting oscillatory terms) |<£>o(R)| 2 ~ {SEp" ' wn i° n 
implies that 

]T MR)| 2 iV 2 JdR-±- ~ 27riV 2 InL, (14) 

diverges as the linear size of the system L — > +00. If we insist in normalizing the zero-mode wavefunction to unity, 
the normalization constant Nq ~ (27rlni) _1 ' /2 — > as L — >• +00. This implies that the hybridization of the zero-mode 
with the localized a level scales as ~ (InL) -1 / 2 . Though this vanishes in the thermodynamic limit, it does so very 
slowly. For instance, for a typical experimental device ao x L 1000 nm. Since ciq ~ 0.25 nm, L ~ 4000,which implies 
that iV 2 ~ 1/50 = 0.02 (i.e. N ~ 0.1). 




n 

FIG. 5: sp 2 -sp 3 hybridization of the dangling a orbital in the apical atom (cf. Fig. [I] in the main text) with the it orbitals on 
the same atom vs. angle between the dangling orbital and the Graphene plane. 
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At this point, it is worth considering also the interaction terms between the electrons in the dangling a level and 
the electrons in the 7r band. Those terms result from the Coulomb interaction, 

where p(r) — v I'J l (r)\E ,Q (r). If we restrict ourselves to the subspace spanned by the a level and the ir and zero mode, 
then we can approximate: 

* Q (r)«^(r)d Q +^ (r)<. (16) 

Hence, 

p(r) = Mr)| 2 4d Q + Mr)|M < + ^(r)^ (r)4< + ^(r)^(r)vrt d« (17) 



Upon substitution of the above in Eq. (151, and upon neglecting terms that involve inter-level transitions, we obtain: 

flint = U<j7rn d n Q + U^nl - J H d^TT^TTg d 13 , (18) 
where we have subtracted the self- interaction (i.e. Hubbard-?/) term for the a level (~ Un^). 

U ™=^J dr|^(r)| 2 MrO| 2 ^-L^, (19) 

U mr =- f dr Mr)| 2 Mr')| 2 ^^, (20) 
e J |r- r'| 

J H = e -J dr^(r)^(r)^o(r')^(r')^7| ■ ( 21 ) 

We also recall that, on dimensional grounds, 

e 2 

U ~ — (22) 

Before trying to estimate the above Coulomb integrals, we note that, by using the following identity for the Pauli 
matrices 

E (0?(*°)S = 2a?$-W. ( 23 ) 

a— x,y,z 

we can recast 

«£*{o*o <*" = 2S <* • S " + ^ n "> (2 4 ) 



where S% = 7rJ Q tt£ and S v = 4a Thus, 



fldTr = Udnnd.no + U^nl - JffS d • S , (25) 

where 

t^r =tJdn- J -f (26) 

Let us next find the scaling of the different couplings U^, U ad) and J// with the system size. To this end, we recall 
that (see e.g. Ref . H3f : 

po(r) = iVo^ cos((K + K') • r) sin((K - K') • r + 9 T ), (27) 

e — A d |r| 

Mr)=N d ——, (28) 
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d 



We shall not perform the integrals in Eqs. ( 19|20|2l" ) explicitly but merely note that, from the 



where A, 

dependence on (fo(v) (i.e. on Nq) and, on dimensional grounds, 



IK 



' H 



UdTT 



e 2 


( 1 ) 


ea 


\2irlaL J 


e 2 


1 


ea 


2tt lni' 


2 

e 


1 


ea 


27rlnL' 



(29) 
(30) 
(31) 



Note that scales faster to zero as L — > +oo than the other couplings because the zero-mode self-interaction is 



strongly suppressed by the fact that its wave-function is only marginally localized. Taking 
we obtain: 

U^jr ~ 1 meV, 
U an ~ 100 meV, 
J„ ~ 100 meV, 
U ~ 10 eV. 



10 cV, and L 1 



10 7 , 

(32) 
(33) 
(34) 
(35) 



However, the hybridization energy of the dangling a orbital with the zero mode is given by VNq. For an angle 
9 ~ 2 degrees, \V\ ~ 200 meV, and N = 0.1 for L ~ 4000, |^|7V ~ 20 meV. Indeed, whereas the Hund interaction 
— JnSd ■ So tends to align the spins of the a level and the 7r zero mode ferromagnetically, the hybridization VNq favors 
an anti-ferromagnetic alignment, but since Jh ^> V 2 /\ed + U a7r \ ~ 1 meV, we can safely assume that the triplet (i.e. 
the ferromagnetism) is favored provided the each level is singly occupied, i.e. rid + n = 2. On the other hand, the 
interaction term cx t/ ffcr tends to favor the emptiness of one of the levels, i.e. 114 + no = 1, which implies that a spin 
doublet is favored. In the following section, we analyze the different possibilities for the total spin of the vacancy. 

The above estimates account for the contribution from the direct Coulomb interaction between the levels localized 
at the vacancy. However, other contributions may also exist. For example, the hybridization with the conduction 
band of both levels resulting from a next-nearest neighbor hopping t' may also lead to an additional (RKKY-type) 
contribution to the Hund's rule coupling Jjj 9J. Furthermore, a large enough t' also shifts energy and broadens the 
zero-mode, which lead to its disappearance into the continuum of scattering states. For the above reasons, in the 
main text we have explored the possible magnetic phases treating Jh and U a7r as unknown parameters. 



THE 'ATOMIC LIMIT OF THE VACANCY 



Let us momentarily set the hybridization with the rest of the it band to zero and consider the spectrum of the 
coupled a and 7r (zero-mode) levels. The latter are described by the following Hamiltonian: 



tine = e n d + e d rid + Un d ^n d ^ + U a7! n d n - J H S d ■ S + V^cL^ttq + V^ir^d" 
where we have introduced 

1 



= VVo(O) ~ VN ~ V 



2ir lni 



1/2 



(36) 



(37) 



In what follows, we are going to determine the ground state of the two level 'atom' described by the Hamiltonian 

\2 



in Eq. (36). To this end, we shall rely on the fact that both n-r = n + rid and = (S + S^) = St(St + 1) and 



— + Sg are good quantum numbers. 



1. Doublet: For tit = rid + n = 1, St = ^, = ±| ; the Hilbert (sub)space of the system is spanned by: 

|i,^ = ±i> = 4|0), 

|2,^ = ±i)=7rt|0), 



(38) 

(39) 
(40) 
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with a = \ {a =1) for S T = \ (S T = — |). The eigcn energies are e± = £d+£o± V( £d ^ £ ") 4 I V ^I ^ Since 
I e<i — £ o| ^ 2|V^ T |, we find that the lowest energy state in th sector has an energy 



(41) 



E (n T = 1, 5 T = -) ~ e d - , r 

2 |Cd — £q| 



2. Singlet: For tit = n d + no = 2, with SV = 
the following three states: 

|1> = 



0. In this case, the Hilbert space of the system is spanned by 



1 

71 



4 7r o^ 4 7r ot 



|0>, 



|2)=44|0), 
|3)=4 r 4j0). 



In this subspace the Hamiltonian Hy becomes a 3 x 3 matrix: 

e d + e + Uvtt + ^f- 
H a7T (n T = 2, St = 0) = | y/2V^ 



'2V a7I 




2e d + U + ^f- 



V2V ( 


2e 



(42) 

(43) 
(44) 



(45) 



3Jh 

4 



where we have used that S d • So = \ [S T — S d — S§] = ^S^ — f . The diagonalization of the above Hamiltonian 
is not possible explicity, but under the assumptions that e > e d and U ^> \e d \, \e d \, the energy of the lowest 
state in this sector is: 



E (n T = 2, S T = 0, S Z T = 0) ~ e d + s Q - 2\V a 



1 



\e d + U - e - U a 



kd + U a7T - e \ 



(46) 



3. Triplet: For tit = 2, S — 1, this space is spanned by (note that the states in this case have different values of 
S z , since Hy commutes with S z , it automatically becomes diagonal in this subspace): 



|5 T = 1,^ = +1)=44 T |0), 
1 



\St — 1, S T 



o> 



V2 



4^4 + 4^ i°) 



\S T = l,S T = -l)=d\4 i \0). 



(47) 
(48) 
(49) 



(S T = l,m\H av \S T = l,m') = S m . m , \^e d + e + U a 
Thus, the triplet consists of three degenerates states with energy: 



Jh 
4 



E (n T = 2, S T = 1) = e d + e + U aiT - 



Jh 



4. Triply occupied doublet: For tit — 3, S — |, and S z = ±| the Hilbert (subspace) is spanned by: 

|l,a) = rf^jo^lol )' 
|2,a)=447rJ Q |0) 

To leading order, the energy of these states is given by: 

E 1 = (l,a\H an \l,a) = 2e Q + e d + 2£/ CT7r , 
E 2 = (2,a\H^\2,a) = e + 2e d + 2£/ CT7r + U, 



(50) 



(51) 



(52) 
(53) 
(54) 



(55) 
(56) 



As in the double case, these energies are corrected by terms of order \y a -n^ l\E\ — E 2 \ in the next leading order. 
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5. Full vacancy singlet: For ut = 4, St = S T — 0, it is described by the state: 

|Full) = 444 7rJ |0), (57) 

with energy E full = 2e d + 2e + 4C/ ff7r + U. 

Under the assumption that U ~ 10 eV 3> £q, |e<j| S> V^ ff , we must compare the spin doublet energy with the triplet 
energy in order to determine the spin state of the vacancy: 

A£ trip i et _ double t = E a (n T = 2, S T = 1) - E (n T = l,S T = \) = e Q + U a7r - ^ + } V<T ^ , (58) 

I 4 | e d — £ 1 

For example, if we choose the position of the chemical potential at the Dirac point and measure the energies with 
respect to it, /i = e = 0, 

A-Btriplet- doublet — U ar; h -, — ■ — r ~ U a7r — , (59) 

4 |£<i + ^(T7r| 4 



However, 



since 1 ^ ~ 1 meV is much smaller than U aT[ ^ Jh ^ 100 meV, according to the estimates made in Sect, 
since U a7r ~ Jh we cannot decide whether the triplet or the doublet are favored as the ground state of the dangling 
CT-7T zero-mode system. Thus, in the following we shall analyze both possibilities. 



SCHRIEFFER- WOLF-LIKE PROJECTION 



In order to integrate out the charge fluctuations, we shall project the Anderson model introduced above onto the 
subspace where the vacancy dangling crdevel and the n zero-mode are either in the triplet or doublet state. This 
requires that (n<j) ~ 1 and that (no) — 1 for the triplet case and (no) = for the doublet case (see Sec. [] for a 
discussion of the validity of these approximations). To this end, we introduce a projection operator, P, onto the 
subspace where = 1 and either n = 1 (triplet) or n = (doublet); Q = 1 — P projects states onto the orthogonal 
subspace where either rid = or n ( i = 2. The effective Hamiltonian projected onto the P-invariant subspace can be 
obtained from the expression: 

H eS (E) = PHP + PHQ E _ l QHQ QHP (60) 

which we can use to perturbatively powers of V integrate out the fluctuations in the occupancy of the dangling a 
level. However, the results, depend on whether the ground state of the dangling a level 7r zero-mode complex is a 
doublet or a triplet. 



Doublet case 



Let us assume that the ground state of the the dangling a level 7r zero-mode complex is a doublet (i.e. U a7r > Jh/%)- 
In that case P is the projection operator onto the doublet state of the vacancy level system. To carry out the 
calculation, we first note that 

PHP = H C + PH aT[ P + PH hyh P, (61) 
QHQ = H C + QH a7T Q + QH hyh Q, (62) 
PHQ = PH hyp Q + PH a7I Q. (63) 

Upon neglecting retardation effects (i.e. setting E = £ + E (n T = 1, St = f) — Sq + e^, where £ is the ground state 



of the chain) as well constant terms and using the identity ( 23 ) , we arrive at 



h k = h c + K P ^ c « c t d Pp _ JY^ Pc i adad y o p (64) 

= H c + J K S d -s +V cl a c% + > Vo ^, (65) 
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where, to leading order in Vq, 



Jk = 2 1 V^o | 2 
1 



V K 



Vo 



1 1 

+ 77 

(d U + e d 

1 1 

+ 



£d U + e d 



2\V \ 2 



1 1 

e^f + U + e d 
1 1 



M U + e d 



>0, 
> 



and 



(66) 
(67) 

(68) 
(69) 



Triplet case 

We can repeat the same calculation as above by assuming that the ground state of the vacancy is a triplet state. 
In this case, P projects onto the triplet state of the vacancy dangling er and it zero-mode levels. Thus, neglecting 
retardation E = £q + Eq{tit = 2, St = 1) = £o + € <i + £q + U a7r — Acting upon with the operator ffhyb then 
produces two kinds of excited states where either n d = (having energy ~ £ + e to leading order in V) or n d = 2 
(having energy ~ £ + 2e d + U + e + to leading order in V). Thus, 



H K = H C - J H P t S d ■ S P t + 2J K P t S d ■ s P t + V cl a c% 



(70) 



In the above expression, P t is the projection operator onto the triplet state. This results holds for temperatures below 



the double-triplet energy splitting 



U a7T . As a consequence of rotational invariance, 



1 



Pt^dPt — 2^ T ' 



(71) 



where St = + So and S T — St(St + 1) being St — 1< Thus, we can write Eq. (70) as a S = 1-impurity Kondo 
Hamiltonian: 



where (recall that U„ 



Jh 



H K =H C + J K S T ■ s + VKclodS, 
< for the triplet to be the ground state): 



Jk = \Vo\ 



\e d + U a ^-^\ U + e d + U an - J -f 



>0, 



\e d + U an -lf\ U + e d + U a « - 3 -f 



> 0. 



(72) 

(73) 
(74) 



In the above derivation, we have used that PtS d ■ SoP* = \Pt [S T — Sq — S^] P t — \ and dropped the resulting 
constant term. 



SOLUTION OF THE SCALING EQUATIONS FOR THE KONDO COUPLING 

Let us next consider the solution of the scaling equations for Kondo coupling: 

dj K dlnA(e) 



d In D d In e 

To find the solution to this equation, let us use the ansatz: 



Jk - J K - 



d\nD d\nD' 



d\nD 



(75) 



(76) 
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and require that our choice of f(D) cancels the first term on the right hand-side of (75), which leads to 

d ( f(D) \ 

I » - 0, (77) 



dlnD \A(e = D) 

and therefore f(D) — CA(e = D), where C is a constant. Upon choosing (7 = 1, without loss of generality, we find 
that g(D) obeys the following equation: 

Mm. m 

which can be integrated from D = Dq to D = T to yield (g(Q) = g(Do)): 



1 1 f T A(e) 



9(T) ff(O) J Do e 
Next we recall that j(T) = A(e = T)g(T) and therefore, 



(79) 



7CT) - (80) 



Hence, 



which, when solved for j(T) yields: 



A(T) _ A(£> ) _ /- J A(e 
J(0) 



(81) 



Do 



A(T)j(0)/A(Dq) 
i _ j(Q) r T 

1 A(D ) Jd q 



Again, we can define the Kondo temperature as the value of T where j(T) = C ~ 1, which form the previous 
expression leads to: 

T "^M-CA(T,,^c< ±. (83) 
D e j(0) Jk 

Introducing A(e) ~ [|e| ln 2 (e/£> )] [H] in the above equation leads to to the condition that Tk In (Tk/Dq) ~ Jk 
quoted in the main text. 



SLAVE-BOSON MEAN-FIELD SOLUTION 



In the present study we employ the slave bosons technique combined with the large -TV mean field treatment of 
the resulting action. The slave bosons technique aims at enlarging the Hilbert space of the impurity by introducing 
a slave boson creation operator and a new fermion f\ that act on a new vacuum state |vac) such that the physical 
empty and singly occupied physical states of the impurity orbital are represented as 

|0)=6 t |vac), (84) 
= /t | vac). (85) 

In terms of the slave boson and the new fermion the original impurity creation operator reads S a = f^b, which 
satisfies canonical anti-commutation relations provided the Hilbert space is restricted to the physical space. The spin 
degeneracy is generalized from SU(2) to SU(iV), with u = 1, . . . N, and we shall consider the U — oo limit with a 
constraint of no double occupancy (for N = 2) generalized to 



fet6 + ^/t/ CT = y (86) 

(T 
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for arbitrary N. The resulting effective action is: 



S cS =f*drf^2 fl{d T + e d )U + tfd r b + iX(r) 



N 



N 



f dr f dT'J2ft(r)b(T)g(r-T'p(r')f a (T') (87) 
Jo Jo 



where the Lagrange multiplier A is a dynamical variable introduced to impose the constraint. 

In the large N limit, the partition function can be approximated by means of a saddle point evaluation of its 
defining functional integral. At the saddle point, the Lagrange multiplier becomes static iX(r) — > (iX) — Xq and the 
boson field is replace by its average b — > (b) = \J N/2ro and tf> —> (tf) = \/Wj2rQ. The values of ro and Ao at the 
saddle point are determined by minimizing of the free energy, which leads to the following mean field equations 



2|V| S 



where Gf is the fermion propagator, which reads: 

Gf(iu n ) = (fl(iuj n )f a (iuj n )) 



1 



(89) 



(90) 



iu n - e d - A - r%V 2 g(iu) n ) ' 

We fix the Kondo temperature TV at the onset of condensation of the slave boson, r = 0. Gf reduces to the standard 



Green's function for a free fermion, and therefore Eq. (|88| leads to 

If 



£d + A 



T 



where = l/(e u + 1) is the Fermi factor. Its solution is given by Ao = — Ed- Thus, Eq. (89) becomes: 



(3 



Ed 



(91) 



(92) 



which we have solved numerically by first computing the local Green's function a the tight-binding model of the 
graphene including nearest neigbor t and next-nearest neighbor hopping t' . 
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